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ABSTRACT 



The analysis of quickly damped transverse oscillations of solar coronal loops using magneto-hydrodynamic seismology allow us to 
infer physical parameters that are difficult to measure otherwise. Under the assumption that such damped oscillations are due to the 
resonant conversion of global modes into Alfven oscillations of the tube surface, we carry out a global seismological analysis of a 
large set of coronal loops. A Bayesian hierarchical method is used to obtain distributions for coronal loop physical parameters by 
means of a global analysis of a large number of observations. The resulting distributions summarise global information and constitute 
data-favoured information that can be used for the inversion of individual events. The results strongly suggest that internal Alfven 
travel times along the loop are larger than 100 s and smaller than 540 s with 95% probability. Likewise, the density contrast between 
the loop interior and the surrounding is larger than 2.3 and below 6.9 with 95% probability. 
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1. Introduction 

The discovery of quickly damped transverse oscillations of solar 



coronal loops was first reported by Aschwanden et al. ( 1999 1 and 
[Nakariakov et al.| ( |1999j ) using Transition Region and Coronal 
Explorer (TRACE) observations. The phenomenon was inter- 
preted in terms of the standing linear magnetohydrodynamic 
(MHD) kink mode of a magnetic flux tube in its fundamental 
harmonic. The cause of the quick damping of the oscillations 
has been attributed to the resonant conversion of global motions 
into localised Alfven oscillations at the tube boundary, because 



of the transverse inhomogeneity of the medium ( Goossens et al. 
|2002| l. The essence of the kink mode has been found to be of 



mixed nature Goossens et al.|( 2009| l with a dominant Alfvenic 
character (IGoosse ns et al.| |2012a). Recent reviews on theoret- 
ical aspects of MHD kink waves can be found in Ruderman & 
Erdelyi (2009jl and Goossens et al. ( 2011| l. Observational proper- 
ties of transverse coronal loop oscillations are presented and dis- 
cussed by |Schrijver et al. (|2002l) and 
MHD seismology (iUchida|fl97d 



Aschwanden et al.|(|2002[ ). 
IRoberts et al ('1984) uses 



inversion techniques to infer difficult to measure physical pa- 
rameters combining theory and observations of MHD waves. 
Coronal seismology applications using transverse loop oscilla- 
tion have been successful in determining parameters such as the 



magnetic field strength (Nakariakov & Ofman 



speed ( Arregui et al. 2007; ,Goossens et al. 



2001), the Al fven 
2008f|2012b| l, the 



transversal density structuring (Verwichte et al. 2006 1, or the 
coronal density scale heig ht ( [Andries et al.j |2005[ [Verth et al.j 
^008] [Arregui et al.||20T3) . 



In a previous paper ( [Arregui & Asensio Ramos] [201 l| l, 
we pursued the Bayesian analysis of individual coronal loops 
with the aim of infeiTing their fundamental parameters. In the 
Bayesian framework the inference is given by a distribution, the 
so-called posterior probability distribution, that is a combination 
of how well the observed data are predicted by the model, the 
likelihood, and our state of knowledge on the unknowns before 



considering the data, given by the priors. In that analysis, we 
demonstrated that the inferred values of the Alfven travel time 
are robust. Additionally, a Bayesian analysis is able to give some 
information (at least put some constraints) on the transverse in- 
homogeneity. A result of the work was the fact that the density 
contrast between the coronal loop and the ambient medium is 
poorly constrained by the observations. The posterior distribu- 
tion changes when different prior distributions are used for this 
parameter. However, if an independent measurement of the den- 
sity contrast is available, we demonstrated that the three parame- 
ters can be accurately inferred from the period and damping time 
of the coronal loop oscillations. 

These previous studies have focused on the inversion of 
physical parameters using measured wave properties for partic- 
ular events, on a one-by-one basis, thus obtaining estimates for 
the local properties of the plasma for each particular event. Also, 
in most of the studies, the wave properties that are used as input 
in the inversions consist of quantities that are obtained upon ma- 
nipulation of parameters obtained at the primary stage of the data 
analysis. Some examples are the use of the period and damping 
of the oscillations that are obtained after a fitting of the measured 
time evolution of the displacement in a sequence of imaging ob- 
servations, or the phase speed of propagating waves derived from 
time-distance diagrams for propagating waves. 

This work represents a substantial step forward along two 
lines. First, we go here to a more fundamental level than "Arregui 
& Asensio Ramos ( 201 1) and use the displace ment curves them- 
selves, as measured by [Aschwanden et al. p002[ ), instead of 
period and damping times. The reason is that the assumption 
of a gaussian likelihood function is more appropriate for the 
displacement curves than for the derived quantities. The noise 
statistics in the derived quantities is very complicated to ob- 
tain, given the non-trivial manipulations that are needed to obtain 
them. 
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Second, a key issue in Bayesian parameter inference is the 
use of prior information that accounts for our state of knowledge 
on the unknowns, before considering the data. This prior knowl- 
edge is usually constructed on the basis of informed guesses 
about, e.g., the values and ranges of variation of physical pa- 
rameters, physical constraints imposed by the model, etc. In this 
work, we compute these priors using observational information 
and obtained from the global analysis of a number of observed 
events. This is done by performing a fully consistent analysis 
of a large number of observations using a hierarchical Bayesian 
framework. In the same way as directly measured properties of 
transverse loop oscillations can be summarised by performing 
histograms, from which quantities such as the mean, median 
or standard deviation can be obtained, the Bayesian hierarchi- 
cal framework enable us to obtain similar information for the 
physical parameters that cannot be directly measured and need 
thus to be inferred. As a result, data-favoured distributions for 
the unknown parameters are obtained. They can then be used to 
construct priors based on our current observational information 
of coronal loops. 

The layout of the paper is as follows. Section 2 presents our 
inference approach based on a Bayesian hierarchical model and 
how this is applied to the observations. The sampling of the pos- 
terior and the marginalization is also discussed in this section. 
Section 3 presents the final results of the paper and we end with 
the conclusions in Section 4. 

2. Hierarchical modeling of coronal loops 

If one has direct observational access to a given physical quan- 
tity (for instance, brightness, displacement, etc.), obtaining an 
estimation of the probability density for this quantity is easily 
achieved by just counting events in bins. If observational uncer- 
tainties can shift events from one bin to another, it is possible to 
use Bayesian schemes to take this and other effects into account 
(e.g., the extreme deconvolution technique of |Bovy et al.||20lT| 
and references therein). However, when one is interested in a 
quantity that cannot be directly observed but has to be inferred 
from observations, the situation is not so straightforward. This is 
exactly the problem we have in our case, because we do not have 
direct access to the physical parameters of the oscillating coro- 
nal loops. We propose a Bayesian hierarchical scheme to solve 
this problem. It can be considered as an efficient way to estimate 
a probability density of an unobserved quantity, obtained from 
many observations of quantities that are non-linearly related to 
the one of interest. 

In summary, in this work we impose parametric shapes for 
the priors for all the parameters of interest and learn the value 
of these parameters from a large set of observed coronal loop 
oscillations. The ensuing final priors with their shape inferred 
from the data summarize all the global information we currently 
have for the physical properties of coronal loop oscillations. 

In the following, we first describe the model used to explain 
the observations and how it depends on the physical parameters 
of interest. Afterwards, we describe the hierarchical probability 
model used to explain the complete set of observations and we 
define the hierarchical priors used in this work. We also explain 
how to efficiently sample the high-dimensional posterior proba- 
bility distribution function. 

2.1. Coronal loops oscillation generative model 

Observing oscillations in coronal loops with the aim of carry- 
ing out magneto-hydrodynamical seismology is a very difficult 



task. After an arduous process that requires a detailed analysis 
of the time evolution of images obtained in coronal lines from 
space missions, the time variation of the displacement, d(t), that 
describes the motion of the coronal loop apex at different time 
steps, is obtained (e.g., Aschwanden et al. 2002[l. In order 



to 



extract information from the time evolution of the displacement, 
this quantity is modeled as a combination of a systematic motion 
of the entire loop and a real oscillatory component. Therefore, 
the generative modelHfor the observations is then: 



d(t) = £/trend(0 + do,c(t) + e(0 + b(t). 



(1) 



where e(f) represents the uncertainty of the amplitude mea- 
surement, while b{t) takes into account the presence of any 
remaining uncertainty produced by non-modeled effects like 
background loops, wrong estimation of the noise variance, etc. 
Concerning the standard uncertainty, we assume it has Gaussian 
statistics, with zero mean and time-independent variance o"^. 
Such a simplification means that measurements at different times 
are completely uncorrelated. Additionally, we use an estimation 
of o-„ obtained directly from the observations. The background 
component is assumed, for simplicity, to be also Gaussian with 
zero mean and time-independent variance crf^. 

The oscillatory component is modeled in terms of a sinu- 
soidal with an exponential decay as follows: 



doscit) = A sin 



27r, 

— (t - to) - cf>o 



exp 



t-to 



Td 



(2) 



where A is the amplitude of the oscillatory part, P is its period, to 
is a reference initial time that is fixed from the observations, 0o 
is the initial phase of the oscillation and t^ is the damping time 
scale. The detrending of the oscillatory motion is very difficult 
to carry out and might depend on an undetermined (potentially 
large) number of factors. For this reason, it is customary to use a 
simple polynomial function, that absorbs all the unknown effects 
( [Aschwanden et al.| [2002) 1: 



'^trend(0 = 2_^">^'~'oy, 



(3) 



1=0 



where the coefficients a, are obtained for each coronal loop and 
the order A^ is adapted to the needed complexity. For simplic- 
ity, we fix the values of the a, coefficients to those obtained by 
jAschwanden et al.]p002| l because no physical information is ex- 
tracted from them. 

The generative model that we have written does not allow 
us to extract much physical information. The period and damp- 
ing time are purely observational quantities and we need to re- 
late them with the physical conditions in the loops. To this end, 
we propose the resonantly damped MHD kink mode interpreta- 
tion of quickly damped transverse oscillations of coronal loops 
( Ruderman & Roberts 2002[ Goossens et al. 2002] ) to explain 
the observed period P and damping time t^. This approxima- 
tion applies to a straight cylindrically symmetric magnetic flux 
tube with a uniform magnetic field pointing along the axis of 
the tube. Under the zero plasma-/? approximation, coronal loops 
can be considered to be density enhancements with a constant 
internal density, p,, a constant external density, p^ < p,, and a 
non-uniform transitional layer of thickness / that connects both 
regions. Following |Goossens et al.| ( [2008j , it is possible to give 



' A generative model defines a parametric description of the signal, 
taking into account the presence of observational uncertainties and its 
statistical properties. 
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the following analytical expression for P and r^ under the thin 
tube and thin boundary approximations: 



ta 



V2 



^+1 



1/2 



and 



P 



!£, + 1 1 

li^- iJJr' 



(4) 



From these considerations, the parameters in which we are inter- 
ested are the internal Alfven travel time, ta, the density contrast 
between the tube and the environment, ^ = Pi/Pe, and the trans- 
verse inhomogeneity length scale in units of the radius of the 
loop, l/R. 

2.2. Hierarchy 

According to the previous model, the oscillatory displacement 
of the /-th coronal loop is determined by the set of parameters 
0i = {ta,^, l/R, a, (po, cTi,}, where we use the vector 0i to compact 
the notation. The Bayesian analysis performed by Arregui & 



Asensio Ramos (2011 1 demonstrated that the constraining power 



of the observations is very limited. Although ta can be success- 
fully estimated from the observations (although with relatively 
large and asymmetric error bars), the situation is much worse for 
the density contrast and the length scale, with the density con- 
trast being the poorest constrained. Arregui & Asensio Ramos 



(201111 have shown that the marginal posterior distribution for 
f is very close to the assumed prior distribution, meaning that 
there is almost no information in the observations to constrain ^. 
The reason why, even in the absence of information for the den- 
sity contrast, the Alfven travel time can be correctly recovered 
has to be found on the specific shapes of the curves in the three- 
dimensional space (ta,^, l/R) pertaining to constant values of P 
and Td, as explained in |Arregui et al.| ( [2007[ ). 

Consider = {0i , ^2, ■ ■ ■ , ^n) to be a vector of length 6A^ that 
contains all the model parameters for all the observed A^ loops. In 
a standard Bayesian approach, the posterior distribution (which 
encodes the updated information about the model parameters) is 
given by: 



P(0|D) 



p(D|Q)p(0) 
P(D) ' 



(5) 



where D - {Di , D2, . . . , D„} refers to the observed data, the mea- 
sured time variation of the displacement, d°^^{t), for all the loops. 
The function p(D\0) is the likelihood, that measures the prob- 
ability of getting a set of observed displacements for a given 
combination of the parameters. Viewed as a function of the pa- 
rameters 0, the likelihood measures the quality of the parametric 
model to explain the observations. Finally, the function p(&) is 
the prior distribution that encodes a-priori information about the 
model parameters, while p(D) is the evidence. Given that p(D) 
does not depend on the model parameters, it is just a multiplica- 
tive constant and can be dropped from the calculations. 

The quantities with physical interest in our problem are ta, ^ 
and l/R. They are obviously directly unobservable. For this rea- 
son, one cannot use the standard histogram to a set of observed 
coronal loops with the aim of obtaining their general physical 
properties. It is widely known that Bayesian hierarchical mod- 
els constitute a very powerful way to overcome this difficulty 
(REF). The idea behind hierarchical models is extremely sim- 
ple. The priors p{Q) used in Eq. (B) are made dependent on a set 
of hyperparameters Q, which are then included in the inference 
scheme. Formally, the posterior is given by: 



pi&,n\D) 



p{Dmp(&\£i)p(n) 
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Fig. 1. Graphical model representing the hierarchical Bayesian 
scheme that we used to analyze the set of coronal loop oscil- 
lations. Open circles represent random variables (note that both 
model parameters and observations are considered as random 
variables), while the grey circle represents a measured quantity. 
The frame labeled "Coronal loop /" represents that the model 
has to be repeated for all the observations. An arrow between 
two nodes illustrates dependency. The nodes that are outside the 
frame are the hyperparameters of the model and are common to 
all coronal loops. 



where we have used the general fact that p(&,£i) = 
p{&\fl)p{£l). Note that we have dropped the dependence of the 
likelihood on Q, given that Q are just hyperparameters or, in 
other words, parameters of the priors. 

If we make the assumption that there is not any correlation 
between any two coronal loops from the set of A^ observations, 
we can largely simplify Eq. (J6]l. In such a case, the likelihood 
and the priors can be factorized, so that the posterior distribution 
simplifies to read: 



p{&, n\D) 



1 



N 



,^,r\ p(Di\edp(0i\mp(m, 



(7) 



where we have made use of standard probability calculus. 

Since the global properties of the physical properties are gov- 
erned by the priors, our aim is to compute the statistical proper- 
ties of their parameters, fl. Consequently, and although it might 
seem counterintuitive, all the individual physical parameters 
are nuisance parameters for us and have to be integrated out from 
the posterior (e.g., Gregory, ^2005): 



/7(n|D) 



Pirn 

MD) i, 



N 

0/ 



d0ip(Di\ei)p(0i\ii), 



(8) 



where we have made used of the fact that the parameters of one 
loop do not affect those of another loop. It is this integration 
operation the one that propagates information from all individual 
loops simultaneously to the hyperparameters. 

2.3. Likelihood 



According to the characteristics of the noise and background 
^^ components of the generative model displayed in Eq. (|T]i, the 
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Table 1. Priors used in this work 



Parameter 


Prior 


Range 


IIR 


Truncated Gaussian 


[0,2] 


^ 


Shifted inverse Gamma 


[l,oo) 


Ta 


Inverse Gamma 


[0,oo) 


A 


Modified Jeffreys' 


[0,oo) 


h 


Uniform 


[-7r,;r] 


o-b 


Modified Jeffreys' 


[0,cx>) 


a 


Uniform 


[0,oo) 


P 


Uniform 


[0,oo) 


y 


Uniform 


[0,oo) 


6 


Modified Jeffreys' 


[0,oo) 


e 


Uniform 


[0,oo) 


'7 


Modified Jeffreys' 


[0,oo) 



likelihood function is a Gaussian. Given that both e{t) and b{t) 
follow the same Gaussian statistics with zero mean although 
with different (time-independent) variances, the total likelihood 
for an individual coronal loop is given by: 



p(A|0,) = Cexp 






{d{tj) - fti-end(fj) - doscitj)) 



2(crl + (rl) 



(9) 



where m, is the number of time steps measured for the /-th loop 
and 



C^(2K)-"'-'^{crl + crl) 



-mi/2 



(10) 



Of importance is to have a good estimation of cr„, the vari- 
ance of the noise. According to jAschwanden et al.| ( |200 2), the 
process of obtaining the time evolution of the displacement for a 
given coronal loop is indeed quite complicated. For this reason, 
we take a conservative approach and use cr„ equal to 10% of the 
maximum absolute displacement in each coronal loop. Our re- 
sults demonstrate that this number is indeed a lower limit to the 
actual uncertainty. 

2.4. Priors 

In the hierarchical Bayesian scheme, as important as the defi- 
nition of the likelihood is the definition of suitable priors. As 
described in the introduction, the idea is that, since the hyperpa- 
rameters of the priors are learnt from all the data simultaneously, 
the resulting prior distributions will be then adapted to the data. 
As a consequence, the prior distributions defined hierarchically 
are generalizations of the standard calculation of a histogram for 
quantities that cannot be directly observed, like l/R, ta and ^. 

To this end, it is favorable to use general probability distri- 
butions that naturally fulfill the boundaries for all the parame- 
ters. The first step is to consider the range of variation of the 

( |2008) , 



model parameters. After Goossens et al 



the model parameters have to fulfill 



we know that 



llRe{0,2], ta>0, ^>1, <^o e [-7r,7r], A > 0. (11) 

Additionally, ta and ^ have upper boundaries that do not emerge 
from the theory but can be estimated based on physical argu- 
ments. We use T^^ = 1500 s and ^™"'' = 100, although their 
precise values are of reduced impact in the final result provided 
that they are large enough. 

The graphical model describing the hierarchy that we con- 
sider in the analysis of coronal loop oscillations is shown in Fig. 




Fig. 2. Examples of the IG(y, 6) distribution, which is a very gen- 
eral distribution for a positive definite quantity. 

[T] The selected priors, which depend on the set of hyperparam- 
eters SI - {a,/3,y,6,e,ri}, are summarized in Tab. [T] We give 
more details in the following. 

2.4.1. Prior for///? 

The theory says that the transverse inhomogeneity length scale 
has to lie in the interval [0, 2], so it is advisable to use a prior that 
automatically fulfills this restriction. We have used a truncated 
normal distribution, that depends on two parameters, a and /3, 
and is given by: 



TN(l/R;a,/3) 



( ^/2^/3r' exp [-(l/R - afldfi')] 0<l/R< ^^^^ 
otherwise. 



Another option that gives very similar results and also depends 
on two parameters is the scaled Beta prior, defined as 



Betail/R;a,/3) 



B(a,/3) 



(l/R)"-\2-l/Rf-\ 



(13) 



where B{a,/3) is the beta function (e.g.,|Abramowitz & Stegun 



19721, which can be computed in terms of the gamma function 
as: 

na)r(J3) 



B(a,/3) = 



r(a+/3)- 



(14) 



2.4.2. Prior for ta 

The Alfven travel time is defined in the interval [0, oo). A quite 
general distribution that is naturally defined in this interval is the 
inverse gamma distribution, which depends on two parameters, 
-y and 6: 



lGiTA;7,S) ^ J^T-""-' 



r(r)- ^'^P'- 



(15) 



This distribution has the advantage of describing variables with 
skewness with only two parameters. The selection of the inverse 
gamma distribution is somehow arbitrary and other distributions 
like the gamma distribution can be chosen. We have verified with 
a few of them that the results are very robust to the precise se- 
lection of the functional form, provided they have sufficient gen- 
erality. A few examples of the shape of this prior are shown in 
Fig. 11 
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2.4.3. Prior for ^ 

The density contrast is a parameter defined in the interval [1, co) 
and scarce information is available as to what the upper limit can 
be. For this reason, we choose a shifted inverse gamma distribu- 
tion, defined as 



SIG(^;e,/7) = ^(^- ir-'exp --^ 



2.4.4. Prior for o-/, 



(16) 



The standard deviation of the background contribution, cr;,, is 
inferred from the data. Given that it is a scale parameter, it is 
customary to use a Jeffreys' prior. Given that o-j, is defined in the 
interval [0, oo) and the Jeffreys' prior is not proper and not well 
defined at zero, we propose a modified Jeffreys' prior ( [Gregory] 
|2005] l: 



MJicr.-^rl^rT)- 



(o-fc + cr°)ln 



o-o + o-™^" 



(17) 



This prior behaves as a Jeffreys' prior (i.e., as crj^^) for cri, » cr'^^ 
and as a uniform prior for cTb <s cr". Consequently, the transi- 
tion parameter cr° is a lower boundary of the Jeffreys' prior. We 
choose the small value cr° = 0.1. We made sure that this value is 

h 

sufficiently small so that the posterior for this parameter peaks at 
larger values and is therefore not influenced by its actual value. 
Concerning cr™", it is made to be very large and its influence on 
the final results is negligible. 



2.4.5. Prior for 0o and A 

Without any additional a-priori information, we choose a flat 
prior for the phase of the oscillation in the interval [-n, n]. This 
uniform prior equals (2?:)"' if -tt < ^o ^ ^ and zero elsewhere. 
The amplitude of the oscillation is a scale parameter that is de- 
fined in the interval [0, oo). For this reason, we choose a modified 
Jeffreys' prior: 



MJ{A-A",A"''''') = 



(a H-A°) In 



A" -I- A" 

A0~ 



(18) 



with A° = 10 ^ (much smaller than the actual amplitude of the 
oscillation) and a very large A™^". 

2.4.6. Priors for liyperparameters 

The hyperpiiors for the hyperparameters a, /3, y, and e are all flat 
in positive real line. For 6 and rj, given that they can be consid- 
ered to be scale parameters, we choose modified Jeffreys' priors 
with very small transition parameter. However, the final results 
are very robust and do not depend on the specific hyperpriors. 



2.5. Sampling the posterior 

It is clear that the integrals of Eq. (jHll cannot be computed analyt- 
ically. Therefore, it is necessary to rely on numerical techniques. 
We cany out the integral using a technique based on a Markov 



Chain Monte Carlo (MCMC; MetropolisetaL 1953 Neal 



1993 I. Instead of the general Metropohs-Hastings method, we 
used a Metropolis-within-Gibbs method ( |Tierney| |1994[ ), which 



has recently been applied by |Sale| ( |2012| l for mapping the ex- 
tinction in the Milky Way using a hierarchical Bayesian modePl 
The reason for using this scheme is that, in principle, the sam- 
pling of the posterior distribution function for every coronal loop 
is independent of the others, except for the presence of the hy- 
perparameters. Therefore, every step of the posterior sampling 
for each coronal loop can be done independently. After one it- 
eration of each chain is carried out, the hyperparameters can be 
updated using a standard Metropolis-Hastings rule. This update 
is then propagated to every coronal loop. The total length of the 
converged Markov chains is of the order of a few hundred thou- 
sands samples. We verified that the Markov chains are converged 
using standard criteria. Finally, the initial 30% of the chain is dis- 
carded to minimize the sample correlation. As well, we use only 
one sample every three to further reduce the correlation. 

2.6. Selection of observations 

Because of the difficulty of observing oscillations in coronal 
loops, some of the cui-ves analyzed by jAschwanden et al.| ( |2002} 
do not really display the behavior that we assume in Sec. |2.1| 
This poses a problem if the generative model of Eq. ([T]) does 
not include the term b(t) because no combination of the model 
parameters yields a fit to the observations whose residual is 
Gaussian with zero mean and variance cr^. However, the inclu- 
sion of (Ti, into the inference solves this issue. The observed 
loops for which the observation is far from a damped sinusoidal 
will display a larger cr/,. 

The total number of coronal loops observed by|Aschwanden| 
|et al.| ( |2002] l is 30. The number of random variables is then 
6A^ + 6 = 186, the model parameters for each loop, including 
the standard deviation of the background contribution, plus the 
hyperparameters . 



3. Results 

3. 1 . Inference about model parameters 

The output of the MCMC code are samples of the model pa- 
rameters which are distributed according to the joint poste- 
rior /:>(0, Q|D). To this, we have to add the advantage that the 
Markov chain for a certain parameter is distributed according to 
the marginal posterior distribution of this parameter, so the in- 
tegrals of Eq. (J8]l are automatically obtained. Figure l3] displays 
the marginal posteriors for a sample of 5 among the 30 coronal 
loops that we consider in this work. 

The upper row shows the marginal posterior for the Alfven 
travel time, which are well constrained in all the cases. The 
marginal posteriors display a conspicuous peak, although the 
confidence intervals are clearly asymmetric. This is similar to 
the findings of Arregui & Asensio Ramos" ("20111), although in 
that paper we did not fit the whole time evolution. 

The second and third rows show the marginal posteriors for 
the density contrast and for the length scale. It is clear from Eq. 
(H I that the length scale and the density contrast are intimately 
related. A fixed value of tj/P can be obtained with an infinite 
number of combinations of ^ and l/R. Therefore, it is almost 
impossible to get reliable information for each parameter sep- 
arately unless a strong a-priori information is available for any 



of the two (see Arregui et al. 2007 Arregui & Asensio Ramos 



|2011[). Our results show a very interesting phenomenon that is a 



- IDL and Fortran 2003 codes can be downloaded from 
https : //github . com/aasensio/mcmc. 
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Fig. 3. Posterior distributions for the model parameters of a sample of five coronal loops. They display the state of knowledge for all 
physical parameters of all loops when the observations are taken into account. The inferred Alfven travel time (first row), density 
contrast (second row), length scale (third row), the standard deviation of the background (fourth row), derived oscillation period 
(fifth row) and damping time (sixth row). The last row shows the original oscillation corrected for the trend (black curve) and the 

best fit (red curve). The black error bars are those associated with cr„, while the red error bars are obtained using (o"^ + cr^) 



direct consequence of the hierarchical scheme. The fact that we 
assume that the priors for ^ and l/R have to be the same for all the 
observed coronal loops introduces a large amount of information 
into the inference. This results into very well defined posteriors 



both for the density contrast and the length scale. The strong 
constraint imposed by the hierarchical model induces that the 
density contrast is roughly the same for all loops, and the trans- 
verse inhomogeneity length scale is the one changing from loop 
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Fig. 4. Inferred values for the parameters (hyperparameters) that define the assumed probability distribution functions for IjR, ta 
and ^. The hyperparameters a and /3 define the prior for l/R, j and 6 are used for the prior for ta and e and r\ define the prior for ^. 



to loop. We conclude that, under the assumption that the phys- 
ical properties of all coronal loops are extracted from common 
probability distribution functions, the damping time scale is fun- 
damentally determined by the transverse inhomogeneity length 
scale. 

The fourth row shows the information inferred for the stan- 
dard deviation of the background component. Interestingly, a-\, 
is always non-negligible, meaning that none of the observed 
coronal loops displays a pure damped sinusoidal oscillation. 
Additionally, the distribution is very well defined in all cases, so 
that it is possible to reliably characterize this background com- 
ponent. 

Although Ta, ^ and IjR are the physical parameters behind 
the model, it is possible to compute the marginal posteriors 
for derived parameters. Using Eq. (J4|, we have computed the 
marginal posteriors for the period and the damping time, which 
are shown in the fifth and sixth rows of Fig. [3] An interesting 
property of these posteriors is that, although some of the model 
parameters might not be strongly constrained, P and r^ are very 
well constrained from the observations. The marginal posteriors 
are really close to Gaussian, which reinforces the assumption 



used in Arregui & Asensio Ramos (201 1 1 of a Gaussian likeli- 
hood with diagonal covariance matrix. 

Finally, the lowest row of Fig. l3] displays the measured dis- 
placement for each loop and the best fit (roughly equivalent to 
the least-squares solution, except for the presence of the priors). 
The black error bars are obtained using the estimated value of 
(Jn, while the red error bars are obtained by adding in quadrature 
(Jn and (T[,. 



3.2. Global properties of coronal loops 

The hierarchical structure of our model allows us to obtain the 
general properties of coronal loops. To this end, we show in Fig. 
fflthe inferred distributions for the hyperparameters that describe 
the prior distributions described in 52.4 The first column shows 
the results for a and /3, that are the parameters of the truncated 
normal distribution for l/R. The results indicate that these hy- 
perparameters have very well defined values. The median values 
are amed ~ 0.85 and ySmed ~ 0.36. Likewise, the results for the 
hyperparameters of the prior for ta are also well defined, with 
Tmed ~ 4.4 and ^med ~ 870. The situation is less favourable for 
the hyperparameters of the prior for ^, probably a consequence 
of the fact that a single inverse gamma distribution is not able 
to capture the complexity of the global properties of f over the 
whole sample of coronal loops. 

Once the hyperparameters are known, it is possible to use 
this information to get the global properties of the physical prop- 
erties of coronal loops. The first approach is to follow what is 
known as the type-II maximum likelihood approximation. In this 
case, we simply evaluate the parametric priors defined in 5 2.4 at 



the most probable values of their parameters, obtained from the 
peaks on Fig. |4] The results are shown as blue lines in Fig. 5] 
Another way, that fully takes into account the presence of un- 
certainties in the hyperparameters, is to use the A^, Monte Carlo 
samples of a, /3, y, 6, e and rj from the posterior to evaluate the 
following marginalized distributions: 

(pd/R)) = —yTN(l/R;ai,Pi,0,2) 

^" /=i 

1 '^■' 
ipijA)) = —YlG{TA,yi,dd 
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Fig. 5. Inferred distributions for the transverse inhomogeneity length scale (left panel), the Alfven travel time (central panel) and the 
density contrast between the tube and the environment (right panel). Grey curves represent the marginalized inferred distribution, 
obtained as the mean of the priors of 5 2.4 with parameters distributed according to Fig. [4] Blue lines are the distributions of 5 2.4 
evaluated at the peak of the distributions of Fig. HI 



ipm = -VsiG(^;e,-,7?,). 



(19) 



The previous expressions are the Monte Carlo approximations to 
the marginalization of the hyperparameters from the hyperpriors. 
These distributions are shown as grey lines in Fig.lS] 

The distributions shown in Fig.|5] which constitute the main 
result of this paper, represent the underlying distribution from 
which the values of l/R, ta and ^ have been sampled, under the 
assumption that this global distribution is shared among all the 
coronal loops. Consequently, they are generalized histograms of 
these unobserved quantities, which akeady take into account any 
possible degeneracy and uncertainty during the inference pro- 
cess. They represent a data-favored updated prior for the param- 
eters of the model. These priors can be used in the future when 
making seismological analysis of coronal loops using the reso- 
nantly damped magneto-hydrodynamic kink mode interpretation 
of quickly damped transverse oscillations. 

Concerning the transverse inhomogenity length scale, the 
left panel of Fig. l5] demonstrates that roughly all allowed values 
are possible. However, the slight shift of the distribution shows 
that there is a small preference for l/R < 1. Concerning the 
Alfven travel time, it is clear from the central panel of Fig. l5] 
that the most probable value for ta is ~ 160 s, with a median 
value of ~ 212 s. The Alfven travel time is below ~ 540 s are 
and above ~ 100 s with 95% probability. In a surely oversimpli- 
fied situation in which the typical length L and density p of the 
coronal loop is known with precision, the Alfven travel time lim- 
its that we obtain might be used to put some general constraints 
on the magnetic field. Given that: 



Ta 



L 

va 



B' 



(20) 



where va is the Alfven velocity and B is the magnetic field. For 
instance, if L ~ 100 Mm and p ~ 10"'^ g cm"^, we end up 
with 6 G < B < 35 G. If the density is an order of magnitude 
larger, the magnetic field range increases in a factor VlO. The 
most probable value of the magnetic field, corresponding to the 
peak of the Alfven travel time in Fig.lSltums out to be B ~ 16 G. 
These figures are just plain estimations based on an unrealistic 
situation in which the properties of the coronal loop are known. 
The information gained for the density contrast is also very 
interesting. We remind that the strong constraint for this param- 



eter is a direct consequence of the hierarchical scheme, which 
forced the same distribution for all observed coronal loops. 
According to our results, the density contrast is above 2.3 and 
below 6.9 with 95% probability, with a median value of 3.8. 

Finally, we display in Fig. l6] the comparison between our 
results and what one would obtain using a simple histogram 
with the inferred value of the parameters. To this end, we have 
used the inferred values of ta and l/R that were obtained by 



Arregui & Asensio Ramos (2011), complemented with the re- 
sults of applying the Bayesian formalism presented of |Arregui &| 
[Asensio Ramos ( 2011) to the observations collected in Table 1 
of Verwichte et al. (12013). A Jeffreys' prior in the range [1.2, 50] 



is used for the density contrast and an uncertainty of 10% is used 
for the period and damping time if no measurement is available. 
Although the results are somehow comparable, note that the er- 
ror bars are not taken into account in the histogram. This is of 
special relevance for l/R and less important for ta, where the 
inferred values are less uncertain. 



4. Conclusions 

This paper presented the inference of the global physical prop- 
erties of coronal loops obtained through MHD seismology. We 
have obtained the inferred distribution of the Alfven travel time, 
size of the transition layer between the surroundings and the 
coronal loop density enhancement. These distributions are valid 
under the assumption that the properties of all coronal loops 
are just realizations of some underlying distributions. The re- 
sults demonstrate that sharp transitions between the surround- 
ing and the internal media are slightly favored. Additionally, we 
have found that Alfven travel times are in the interval [100, 540] 
s with 95% probability. If the length and density of the coro- 
nal loop are known, this poses some constraints on the magnetic 
field strength in the loop. Likewise, the density contrast between 
the loop interior and the surrounding is in the interval [2.3,6.9] 
with 95% probability. 

Our contribution improves over our previous approach. First, 
we make the model closer to the observation, by using a gener- 
ative model to explain the measured displacements. Second, we 
use a method that obtains global information for model param- 
eters that cannot be directly measured but need to be inferred. 
Our results allowed us to construct informative priors that can 
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Fig. 6. Comparison between the inferred distributions shown in Fig. p and a simple histogram carried out with the inferred values 
of Tyi and IIR using the formalism of,Arregui & Asensio Ramos, (,201 Ij. 



be used for inversions of individual events. The inference then 
takes into account prior beliefs extracted from data. 

Apart from the extraordinary difficulty of extracting the os- 
cillations in coronal loops, the potential to massively apply 
MHD seismology techniques is now larger than ever thanks 
to the continuous observations of the Atmospheric Imaging 
Assembly (AIA; |Lemen et al.[ |2012| l on board the Solar 
Dynamics Observatory (SDO). Its high-temporal cadence of 12 s 
and high spatial resolution of ~0.6 arcsec make them the perfect 
instrument to follow these oscillatory events and extract reliable 
physical information from these coronal events. 
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